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Abstract. We measure the flow of granular materials inside a quasi-two 
dimensional silo as it drains and compare the data with some existing models. 
The particles inside the silo are imaged and tracked with unprecedented resolution 
in both space and time to obtain their velocity and diffusion properties. The data 
obtained by varying the orifice width and the hopper angle allows us to thoroughly 
test models of gravity driven flows inside these geometries. All of our measured 
velocity profiles are smooth and free of the shock-like discontinuities ("rupture 
zones") predicted by critical state soil mechanics. On the other hand, we find 
that the simple Kinematic Model accurately captures the mean velocity profile 
near the orifice, although it fails to describe the rapid transition to plug flow far 
away from the orifice. The measured diffusion length b, the only free parameter in 
the model, is not constant as usually assumed, but increases with both the height 
above the orifice and the angle of the hopper. We discuss improvements to the 
model to account for the differences. From our data, we also directly measure the 
diffusion of the particles and find it to be significantly less than predicted by the 
Void Model, which provides the classical microscopic derivation of the Kinematic 
Model in terms of diffusing voids in the packing. However, the experimental data 
is consistent with the recently proposed Spot Model, based on a simple mechanism 
for cooperative diffusion. Finally, we discuss the flow rate as a function of the 
orifice width and hopper angles. We find that the flow rate scales with the orifice 
size to the power of 1.5, consistent with dimensional analysis. Interestingly, the 
flow rate increases when the funnel angle is increased. 



PACS numbers: 45.70.Mg, 66.30.-h 



Submitted to: J. Phys.: Condens. Matter special issue on Granular Media 



Velocity profile of granular flows inside silos and hoppers 



2 



1. Introduction 

Granular materials display a surprisingly complex range of properties which make 
them appear solid or liquid like depending on the applied conditions |2| ■ Because 
the interaction between the grains is dissipative and the thermal energy scale is small 
compared with the energy required to move grains, such materials quickly come to 
rest unless external energy is supplied constantly. Although vibro-fluidization and 
tumbling [3] is frequently used to excite granular materials, flows driven purely by 
gravity can occur in nature as well. Typical granular flows are dense and a fundamental 
statistical theory is not available to describe their properties. One reason for this 
situation is the lack of quantitative data which can be used to test and develop models 
of dense granular flow. In this paper, we focus on flows inside silos and hoppers in 
order to elucidate the nature of the flow and to test existing models. Such systems are 
ubiquitous due to the need to store and process granular materials in devices ranging 
from simple hour glasses to sophisticated nuclear pebble reactors @| E] . 

Several aspects of granular drainage have been studied over the years. Beverloo 
thoroughly investigated the relation between the orifice size and the mass flow rate 
in cylindrical silos and proposed a formula describing the observed dependence 
Using radiography, Baxter, et al observed the density wave in the hopper flow and 
showed various patterns of the wave depending on the particle roughness and the 
hopper angle [7]. 

The velocity field of the flow inside a silo has been described by two different 
approaches. One is based on the critical-state theory of soil mechanics which relates 
stress and density to predict velocity field or mass flow rate |H1 E]- Although this 
approach has the appeal of starting from mechanical considerations, some questionable 
assumptions are made to resolve indeterminacy in the stress tensor, and the resulting 
equations are mathematically ill-posed and can lead to violent singularities |10t ilf| . 
The solutions available for hoppers possess shock- like velocity discontinuities ( "rupture 
zones") in], which are not seen in our experiments (see below). 

The second approach ignores the stress field and attempts a purely kinematic 
description of the velocity profile, starting from an empirical constitutive law. A 
theory of this type was first discussed by Litwiniszyn, who introduced a stochastic 
model in which particles perform random walks through available "cages" E| . 

Later, MuUins independently proposed an equivalent stochastic model of the flow in 
terms of "voids" and extensively developed the continuum limit, where a diffusion 
equation arises Decades later, Caram and Hong revisited the Void Model and 

implemented it explicitly in computer simulations on a triangular lattice (where the 
voids are simply crystal vacancies) '17'. 

As an alternative to the microscopic void picture, Nedderman and Tiiziin derived 
the same continuum equation starting from a constitutive law relating horizontal 
velocity and downward velocity gradient ^1 E| ■ Regardless of its derivation, the 
Kinematic Model predicts velocity fields with only one free parameter. In light of 
its simplicity, early experiments on silo drainage were viewed as successes of the 
model 121 1201 I21| . even though it has since fallen from favor in engineering 
Although the free parameter has been observed to be proportional to grain diameter 
in all experiments, the constant of proportionality does not agree ^| 1211 I22j . 
Furthermore, Medina, et al 23 have reported that the kinematic parameter varies 
within a silo when the flow is analyzed in detail by particle image velocimetry. 

In addition to the studies of the flow pattern, the diffusion of particles has been 
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investigated as well. Hsiau and Hunt and Natarajan, et al 25^ imaged tracer 
particles in a dense flow inside a vertical channel with various boundary wall condition 
to investigate the concept of "granular temperature". From an analysis of velocity 
fluctuations, they found that particles shows normal diffusion and that the diffusivity 
in the stream wise direction is higher than in the transverse direction. Later, Menon 
and Durian used diffusing-wave spectroscopy to measure the dynamics of 100/im glass 
beads inside a three-dimensional flow with improved temporal resolution, albeit at 
rather small length scales They reported that the particles show ballistic flight 
between collisions over a short time scale, and normal diffusion over longer time scale, 
although the collision distance of 28 nm (1/10,000 of a grain diameter) could perhaps 
be associated with sliding or rotating asperities in frictional contacts. In any case, the 
randomizing gas-like collisions assumed in kinetic theories [23 EHl OH] have not been 
confirmed in any experiments on dense flows. 

With rapid advances in high-speed digital imaging technology, it is now possible 
to simultaneously record thousands of individual particle positions with high spatial 
and temporal resolution. In a recent experiment by our group using this technique, 
the dynamics of 3 mm glass beads near a transparent wall in a three-dimensional 
silo was observed to be sub-ballistic but super-diffusive over short time intervals, and 
diffusive over long time intervals j3()| . The data was argued to be consistent with slow 
cage rearrangement with particles remaining in long-lasting contacts by showing that 
the diffusion scaled only with distance traveled. 

Therefore inconsistencies can be noted in reported results which need to be 
resolved with thorough investigations. In the next section, we introduce the kinematic 
description of granular flow in silos and hoppers in detail. Since the stress field is not 
measured by imaging techniques, we do not assess critical-state mechanical models, 
aside from seeking the presence of the predicted shocks in the velocity field. Then, 
we outline the experimental setup in section [SJ and compare the prediction of models 
with our experiments in section 01 We discuss the implications of the comparison in 
section Eland finally summarize the results in sectional 

2. Models for the mean velocity profiles 

A simple kinematic description of the mean velocity profile in silos and hoppers has 
been developed by since the 1950s, from a variety of theoretical perspectives The 
continuum Kinematic Model starts from an empirical constitutive law relating velocity 
components JHIj which can be derived as a continuum limit of the (earlier) Void 
Model ^lEj- The latter is a more complete theory, because it provides a microscopic 
mechanism for flow, which can be checked by experiments on diffusion and mixing. 
Recent experiments, however, have firmly rejected the void hypothesis. On the other 
hand, an alternative stochastic description, the Spot Model |31ll32j . which starts from 
a cooperative mechanism for random-packing rearrangements, roughly preserves the 
mean flow profile of the Kinematic Model, with much less diffusion and slow cage 
breaking, consistent with experiments |30|. 
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2.1. The Kinematic Model 

Nedderman and Tiiziin proposed a model based on the following constitutive 

law relating velocity components: 

which states that the horizontal velocity u, is proportional to the horizontal gradient 
(i.e. the shear rate) of the downward velocity v. This assumption is based on the 
fact that particles tend to drift horizontally towards a region of faster downward flow 
as they are likely to find more space to move in that direction. Assuming that the 
density fluctuation is small in dense granular regimes, they combined Eq. with the 
incompressibility condition, 

du dv ^ ^^-j 
dx dz ' 
and obtained an equation for the downward velocity, 

dv d^v 

Eq. (PI has a form of the diffusion equation, where time is replaced by the vertical 
coordinate z. When an "initial condition" is given for v at the bottom of the silo at 
z = 0, the velocity diffuses upward. The boundary condition assumed at the side walls 
of the silo is that the velocity is parallel to the wall. Although the authors did not 
discuss this situation, the condition can be naturally generalized to the case where the 
side walls are not vertical. It is written as 

uux — vuz—Q at {x,z) on the side wall, (4) 

where {n^ ,nz) is the normal vector at the boundary. 

For a semi-infinite quasi-two dimensional system (— oo < x < oo) with a point-like 
orifice at 2 = which acts as a source of velocity, a similarity solution exists: 

where, Q is the flow rate per unit thickness of the silo. We refer to the constant of 
proportionality b in Eq. as the "diffusion length," as it is has units of length. We 
provide microscopic understanding of b in section 12.21 

The Kinematic Model has been tested experimentally, and the parameter, b has 
been measured by various groups. Nedderman and Tiiziin observed b « 2.24(i for 
various particle size [19^ . Experiments by MuUins with monodisperse iron ore particles 
imply bra 2d [20] . Medina, et al used the particle image velocimetry (PIV) technique 
to obtain the velocity field and found that the diffusion length increases from b k, l.'ad 
to 6 ~ 4d as the height increases to fit the field. jSHl- Samadani, et al reported b « 3.5(i 
for monodisperse glass beads using difference imaging to find velocity contours [221 • 
All the groups claimed that the prediction of the Kinematic Model qualitatively agree 
with their experiment. The fact that a single fitting parameter b suffices to reproduce 
the entire flow field should be viewed as a major success of the Kinematic Model. 

In order to test the Kinematic Model more thoroughly, we use numerical methods 
to solve the Kinematic Model subject to the same dimensions used in our experiments. 
For this purpose, we define the stream function, ili{x, z) = v{s, z)ds and solve for 
^{x,v) rather than v{x,z). Formulated in terms of ip, the boundary condition turns 
into a Dirichlet one from a rather complicated one given by Eq. Q. Furthermore, 
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it is more convenient for the hopper geometry with inchned boundaries. If the 
width of the system is given by L{z) and the silo is symmetric about its center [e.g. 
—L{z)/2 < X < L{z)/2], the equation and the boundary condition for is given by 

We numericaUy integrate Eq. (jSJ from z = using the Crank-Nicholson method to 
obtain the prediction of the Kinematic Model. 

Due to its continuum formulation, the Kinematic Model cannot predict grain-level 
diffusion and mixing, so we now turn to statistical kinematic models for the velocity 
profile, which postulate mechanisms for random-packing dynamics. 



2.2. The Void Model 

Since Eq. Q has the form of a diffusion equation, where the vertical distance z plays 
the role of "time" , it is clear that any microscopic justification for the Kinematic 
Model should be based on independent random walks. In fact, this is how the model 
was first derived decades earlier, based on statistical considerations. Although the 
continuum approach is more general, in the sense that it is not tied to any specific 
microscopic mechanism, it lacks a clear physical basis, so it is important to consider 
what kind of microscopic mechanisms might support it. 

Litwiniszyn first suggested the idea that particles are confined to a fixed array 
of hypothetical "cages" as they perform random walk from one available cage to 
another during drainage Then, MuUins ^HIEI independently proposed 

an equivalent model in terms of "voids" rather than particles, which is analogous to 
vacancy diffusion in crystals. In his model, particles move passively downward in 
response to the passage of voids, and the voids take directed random walks upward 
after emerging from the orifice. 

Assuming that voids diffuse by non-interacting random walks, it is straightforward 
to show that in continuum limit, at scales larger than the grain diameter, the 
concentration (or probability density) of voids, py, satisfies the diffusion equation, 

dpv /-X 

Since downward velocity v is proportional to the frequency of the void passage, this 
implies Eq. of the Kinematic Model. However, the equivalence of the two model 
assumes that voids can be superimposed without interaction. 

The void model also gives us an interpretation for the kinematic parameter, b. If 
a void undergoes a random horizontal displacement, Axy, while it climbs up by Azy, 
the parameter b is given by 

Var(Ax^) 
^= 2Az. ' 

which is the characteristic length of the void diffusion. However it is very difficult 
to determine b directly from Eq. ©. Axy and Azy cannot be measured from an 
experiment, nor does any a priori choice produce the measured value of b. MuUins 
also deduced b ^ 2d from the velocity profile for round particles {b ~ d/'i for irregular 
particles) without specifying the value of Axy and Azy. By contrast, Caram and 
Hong ^7] assume a void makes an one-to-one exchange with particles on a regular 
lattice when they later revisited the void model. It is noteworthy that any regular 
lattice of hard sphere packing under predicts b {b ^ d) [31] . 
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The void model faces more serious problems when it is used to predict diffusion 
and mixing, which was not done by its proponents. If a tracer particle is placed in a 
uniform flow driven by voids, the particle does a directed random walk downward with 
precisely the same diffusion length as the voids moving up. Thus particles are easily 
mixed before they drop by a few particle diameters, which goes against our everyday 
experience and experiments (see below). 



2.3. The Spot Model 

To address these contradictions, Bazant et al |31[I32| proposed the Spot Model, which 
starts from a mechanism for cooperative diffusion in a dense random packing. It has 
roughly the same mean flow as in the kinematic model, because it also assumes that 
particles move in response to upward diffusing free volume, but this excess volume is 
carried in extended "spots" of slightly enhanced interstitial volume, not in voids. 
The kinematic parameter, 6, is now set by the diffusion length for spots, 
_ Var(Ax,) 

2Az, ' ^' 
where, Aa;^ and Az^ are spot displacements in x and z directions, respectively. Unlike 
a void which is a vacancy capable of being filled by an entire particle, however, a spot 
carries small fraction of interstitial space spread across an extended region and causes 
all affected particle to move (on average) as a block with the same displacement in 
the opposite direction to the spot. 

Of course, there are more complicated internal rearrangements, which can be 
taken into account to achieve accurate spot-based simulations [g3| . but the simplest 
mathematical model already captures many essential features of dense drainage )31l 
1221 ■ For example, it is easy to see that the spot mechanism greatly reduces the 
diffusion length of particles, compared to the diffusion length of free volume. Suppose 
that a spot carries a total free volume Vs, and causes equal displacements (Accp, Azp), 
among Np particles of volume Vp. The particle displacement can be related to the spot 
displacement (Axp, Azp) by an approximate expression of total volume conservation, 

Vp{/^xp, Azp) = -Vs{^xs, Az,), (10) 

which ignores boundary effects at the edge of the spot. From this relation, we can 
compute the particle diffusion length, 

Var(A2;p) u;2Var(Aa;,) 
o» = — Tn. — — — — i = wOs H 

which is smaller than the spot diffusion length by a factor, w = Vs/NpVp. This can in 
turn be related to the change, A(j), in local volume fraction, (p, caused by the presence 
of the spot, 

.^^^^.^ (12) 
b NpVp (j> ^ ' 

It is well known from simulations and experiments that the volume fraction fluctuates 
on the order of 1% in a dense flow, so the Spot Model thus predicts w — bp/b — 
O(10~^). (In our experiments, the local area fraction of glass beads near the viewing 
wall varies by less than three percent.) The estimate of w is further reduced by noting 
that spots occur in large numbers and overlap, so that each spot contributes only a 
small part of the change in local volume fraction. We will test this prediction in our 
experiments. 
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Figure 1. (color online), (a) A raw image of the glass beads acquired with the 
high-speed camera, and (b) the preprocessed image along with the position of the 
centroid of the identified particle ( X ) . 

3. Experimental procedure 

3.1. Experimental setup 

Our experimental apparatus and procedure is similar to that used in our previous 
report "W. We use black glass beads (d = 3.0 ± 0.1 mm) in a quasi-two-dimensional 
silo with length L = 20.0cm (67d), height H — 90.0cm (300c?), and thickness 
D — 2.5cm (8.3c?). The particles near the front wall of the silo are measured through 
the transparent glass. The slight polydispersity reduces the tendency for hexagonal 
packing to occur near the wall. The thickness of the silo D is large enough that finite- 
size effects are not significant. We obtain similar results for both mean velocity and 
diffusion when we increase D |30j. A distributed filling procedure was used to fill the 
silo with the grains. The orifice is opened and steady state flow is allowed to develop 
before acquiring the images used for determining particle positions. 

We view a rectangular region of 20.0 x 50.0 cm above the orifice with a resolution of 
256 X 1280 pixels. Therefore each particle diameter corresponds to c? = 7.7 pixels. The 
images are acquired at a rate of 125 frames per second. The camera memory allows 
2048 consecutive images to be stored at this resolution and therefore the maximum 
interval over which we can track a particle is about 16.4 s. 

For the funnels in the hopper, plexi-glass wedges are placed on top of the bottom 
plate. The surface property of wedge boundaries is identical to the side walls. We use 
wedges with three different angles, 9 — 30°, 45° and 60°. The orifice size W — 18 mm 
is fixed for the hopper experiments while, it is varied with W = 12, 16 and 20 mm for 
the silo. To gain good statistics, three experiments are conducted for each funnel angle 
and orifice size. We also use data from a wider range of orifice sizes than acquired 
during our previous study 30 in section 

3.2. Particle tracking 

To identify the locations of particles from images, we employ the algorithm proposed 
by Crocker and Weeks [35 ■ In this algorithm, the raw images are preprocessed to 
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Figure 2. (color online), (a) Contour plot of the average downward velocity 
field, t) in a fiat-bottomed silo with an orifice width, W = f6mm. (b) d as a 
function of x at the two heights, zi = 9.1d and Z2 = 29. fd indicated with gray 
dotted lines in (a). The result from the Kinematic Model in the same geometry 
fits best with b = l.3d for the zi profile, and b = 2.3d for the Z2 profile. The 
result from the model for the 22 profile with b fitted at zi (narrow solid curve) is 
also shown. 



reduce the noise and the background. This involves convolving the image with a 
Gaussian filter and then an average filter of roughly d pixels respectively. The particle 
location is then identified with the centroid around the local maximum brightness 
pixel in the modified image. To optimize the particle tracking for our experiment, 
the algorithm was also further customized. Because the glass beads are circular, we 
use circular shaped filter. We also set an intensity cutoff to discard the blur images 
of particles located deep away from the front wall. A sample of an image before and 
after the processing is shown in figure ^ The position of the located particles is also 
superposed. 

After particles are located frame by frame, their trajectories should be retrieved 
by "connecting" their positions in time. We associate a particle in a frame with another 
in the next frame which is within a radius of 0.66 d pixels around the original position. 
This simple method works well avoiding more complicated multiple associations except 
very near the orifice where the particles move more than 0.66 d pixels per frame. The 
particles can be tracked there by using a faster frame rate, but we do not do so here 
since bulk flow, and not orifice dynamics, is the focus of our study. 

4. Analysis of the experimental results 

Comparison of the measured velocity profiles with the Kinematic Model 

We first compare the data from the flat-bottom silo with the model. Figure|2Ia) shows 
the contour plot of the average downward velocity v. The mean velocity is obtained 
by dividing the observation window into square cells of size 1.6d x 1.6d. Then in 
each cell, the average is performed over the displacements of all the particles passing 
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Figure 3. (color online), (a) Contour plot of the average downward velocity 
field, j; in a hopper with angle, 8 = 30° , and W = 18 mm. (b) d as a function of 
X at the two heights, zi = 9. Id and 22 = 29. Id indicated with gray dotted lines 
in (a). The result from the Kinematic Model fits best with b = 2.1d for the 21 
profile and b = 2.6d for the 22 profile. The result from the model for the 22 profile 
with b fitted at 21 (narrow solid curve) is also shown. 
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Figure 4. (color online), (a) Contour plot of the average downward velocity 
field, I) in a hopper with angle, = 45° , and W = 18 mm. (b) w as a function of 
X at the two heights, 21 = 9. Id and 22 = 29. Id indicated with gray dotted lines 
in (a). The result from the Kinematic Model fits best with b = 2.1d for the 21 
profile and b = 2.6d for the 22 profile. The result from the model for the 22 profile 
with b fitted at 21 (narrow solid curve) is also shown. 
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Figure 5. (color online), (a) Contour plot of the average downward velocity 
field, D in a hopper with angle, = 60° , and W = 18 mm. (b) j; as a function of 
X at the two heights, zi = 9. Id and Z2 = 29. Id indicated as gray dotted lines in 
(a). The result from the Kinematic Model fits best with b = 2.6d for the zi profile 
and b = 3.2d for the Z2 profile. The result from the model for the Z2 profile with 
b fitted at zi (narrow solid curve) is also shown. 

through the cell. We again take the average of the field from three experiments. The 
data across experiments shows little variation, which confirms that the velocity field 
is well-defined and stationary. Thus we do not show the error bar in the plots of this 
paper unless the concerned quantity has visible fluctuations. 

The contour plot shows that v is maximum right at the orifice and appears to 
"diffuse" upward, in qualitative agreement with the models discussed above. The 
regions in the left and right corner made by the side walls and the bottom plate remain 
stagnant, and the boundary of mobile region has a parabolic shape. In figure |21^b) , 
we show the profiles v{x) at two cross sections zi — 9. Id and Z2 = 29. Id [dotted lines 
in figure |2Ia)] with the fit to the Kinematic Model. The diffusion length, b — 1.3c? 
was the best fit for the profile at zi. However, 6 becomes larger when z increases as 
some previous reports have also reported [51 . The profile at Z2 is best fit with 
b = 2.3d, but it has a flattened shape at center with thinner tail indicating further 
obvious deviations from the model. 

The velocity profiles from the experiments with different orifice width turn out 
to coincide when they are normalized by the flowrate as is commonly observed in 
other dense granular flows . Thus the best fitting value of b is independent of the 
flowrate. The dependence of the flowrate on the orifice width willl be discussed in the 
next subsection. 

We performed similar analysis of the experiments with the hoppers. The contour 
plots along with the profiles at z — zi, Z2 for the angles, 30°, 45° and 60° are presented 
in figure|21Eland[Slrespectively. As angle is increased, the stagnant region is diminished 
as the particles slip on the wedge. At zi, the critical angle over which slip occurs is 
between 30° and 45° and at Z2, it is between 45° and 60°. However the shape of equi- 
velocity contours well above the funnel is not affected significantly by the funnel's 
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Figure 6. (color online). The dependence of flow rate on (a) the effective orifice 
width, W/d — 1 in a flat-bottomed silo (log- log scale), and (b) the funnel angle 9 
in a hopper with a flxed orifice width. The flow rate is measured averaging the 
downward velocity in the plug-flow region. The fitting of (a) validates the result 
of a dimensional analysis, Q oc (W — d)^'^ . 



detailed shape. 

The value of b to obtain the best fit depends on the angle of the hopper as well. It 
increases from 2. Id to 2.8d for zi, and from 2.6c? to A.5d for Z2 as the angle is increased. 

Although we observe some quantitative discrepancies with the simple Kinematic 
Model with a constant coefficient, b, the flow is at least qualitatively consistent. This 
appears not to be the case with continuum models from critical-state soil mechanics Wl , 
which generally predict sharp, shock-like discontinuities in velocity (and stress, which 
we do not measure) within the silo, especially near corners. We see no such abrupt 
jumps in velocity in the silo, only rather smooth velocity profiles. 

4-2. Flow rate dependence on the orifice size and the funnel angle 

The mass flow rate in a silo during discharge was an important subject of early 
research. Using drainage experiments in cylindrical silos with a circular oriflce, 
Beverloo, et al |H| reported a relation known as the Beverloo correlation 

Q (X py/g {W - kd)^-^, k = lA (13) 

where Q is the mass flow rate, p is the bulk density of packing, g is the gravitational 
constant and W is the diameter of orifice. It is usually argued that Q oc {W—(^f ^^ 
is the only form which can be deduced from the dimensional analysis as {W — d) is the 
effective diameter (or width) where particle centers can be placed within the orifice, 
but arching and other effects could also introduce the particle diameter d and thus 
another dimensionless parameter, d/W . Instead, the Beverloo correlation includes 
a somewhat controversial factor W — kd, where the empirical factor k is claimed to 
derive from the region near the orifice rim which obstructs the passage of particles. 
This picture could be consistent with the concept of an "empty annulus" proposed by 
Brown and Richards |3f)| . 

For a slit orifice with a quasi- two-dimensional silo as in our experiment, the 
dependence can be obtained to be 

Q(xpy^{D-d){W-kd)^-^, (14) 
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because the flow rate is linear with systera depth D — d |. We investigated the flow 
rate dependence on orifice width using our data. Although the discharged mass flux 
is not directly measured, we use the overall average velocity, v* = Q/L to obtain the 
flow rate. 

Figure Ufa) shows the relation between the flow rate and orifice size in log-log 
scale. When = 1, the data fits to a power law scaling with an exponent of 1.48. 
Although k = 0.94 gives the exact exponent of 1.5, we do not attach much importance 
to the deviation as our flow rate measure is indirect. However, it is sufficient to check 
that the Beverloo correlation (dimensional analysis) holds in a 2-D silo. 

We also investigated how the funnel angle affects the flow rate. In order to 
compare the rate at a fixed orifice width, wc interpolate the rate with W — 18 mm 
from data with W = 16, 20 and 24 mm for the silo experiment. Figure El^b) shows 
a consistent increase in the flow rate as the angle increases. The flow rate in the 
60° funnel turns out to be about 33% more than that in the flat-bottom silo. This 
dependence is consistent with the data from Ref. [37], although the reported increase 
of the flow rate is smaller than our data. We believe the increased flow rate largely 
comes from the fact that the smooth rigid boundary facilitates the passage of particles. 
As seen clearly from figureEtb)'-4nib), the stagnant zone present in the corners of flat- 
bottom silo gets replaced by wedges. Thus the slip velocity at the boundary increases 
as the angles increases, which make the out-going flow at the orifice (z — 0) more 
uniform and shear-free. This effect appears to allow the particles to exit the orifice 
more easily. 

4-3. Diffusion of particles in an uniform flow 

As explained in section |21 particle diffusion is a key property to distinguish between 
different possible microscopic mechanisms for dense granular flow. The Void Model 
and the Spot Model predict quite similar mean flow profiles (given by the Kinematic 
Model scales much larger than the grain size), but the former predicts bp/b = 0(1) 
while the latter predicts bp/b = 0(10^^). In this section, we briefly discuss 
measurements of particle diffusion in our experiments, as also previously reported 
in Ref. gOI- 

To measure diffusion, the random component of particle displacement is obtained 
by subtracting the average component: 

Ax = Ax'-uAt and Az = Az' + vAt, (15) 

where Ax' (or Az') is the observed displacement in x (or z) direction. Ax (or Az) 
is the random displacement in the same direction, and At is the time gap between 
two consecutive frames and can be increased by any integer multiple. The observation 
window is in a nearly plug-flow region far from the oriflce, where u is negligible and v 
is almost uniform (and set by varying the orifice width). 

The probability density distributions of Ax and Az are observed to display fat 
tails compared to a Gaussian distribution. The statistics of Az also show an anisotropy 
due to gravitational acceleration and inelastic collisions. When the width of the 
distributions is examined as a function of At, the scaling shows a crossover from 
super-diffusion, (Ax^) oc At^-^ and (Az^) cx At^-^, to diffusion, (Ax^) oc (Az^) oc At. 

X It should be noted that the orifice in our system is entirely open from front to back surface. Thus 
the "empty annulus" argument is difficult to apply in the direction of silo depth. We cannot find the 
exact dependence on D because we fix D = 2.5 cm. 
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Figure 7. (color online). The locally measured diffusion length, fe as a function 
of the normalized velocity, v/v*. 



A significant observation is that the hnes of (Ax^) and (Az^) for different v collapse 
into a single line when they are plotted against the distance dropped, vAt, allowing 
us to characterize the dynamics only by distance moved, independent of the flow rate, 
V. We found that this dynamical crossover occurs after a particle falls roughly by its 
diameter irrespective of the flow rate. 

The fact that the dynamics only depends on geometry strongly suggests that 
advection and diffusion have the same physical source (such as a the passage of a void 
or spot). It also suggests that structural rearrangements with long-lasting contacts 
dominate diffusion in dense granular flows, as opposed to ballistic collisions, which 
are central to the kinetic theory of gases. A direct evidence is that the cage-breaking 
length is estimated to be order lOOd from the rate of the nearest neighbor loss [301. Our 
results suggest that the concept of "granular temperature" based on thermodynamic, 
randomizing collisions is of dubious value in slow, dense granular flows. 

Since the free volume models correctly predict the geometry dominated diffusion, 
we can proceed to evaluate them quantitatively. We compute the Peclet number, the 
dimensionless ratio of advection to diffusion, defined as 

2VdAt d , ^ 2VdAt d , 

FCx — hm — - — ^ = - — and Pe^ — um — — ^ = - — (16) 

are interpreted as the distances (in unit of d) for a particle to fall before it diffuses 
by a diameter in a; or z direction, respectively. The large measured values, Pe^: = 320 
and Pez = 150, indicate that advection dominates diffusion. Since b/d k, 2, we also 
find Wx — bp^x/b ~ 1/600 and bp^^/b ~ 1/300, which is consistent with the simple 
prediction of the Spot Model, w « lO^'^ — 10^^. Of course, the data firmly rejects 
the Void Model, which predicts, u; « 1, and cage breaking at the scale of one particle 
diameter. 
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5. Discussion 

In section WJ\ we observed that the Kinematic Model with a constant parameter b is 
not consistent with the experiments. It was found that b depends on the height and the 
funnel angle. In this section, we investigate the validity of two important assumptions 
of the Kinematic Model, namely the constitutive law and the generalized boundary 
condition for the funnel geometry. 

First, we directly check the constitutive law (^) using the results from our 
experiments. In each cell that was used for averaging the velocity, we measure 
the horizontal velocity u, the downward velocity gradient dv/dx, and therefore the 
diffusion length b. Figure [T] shows the distribution of the locally measured values of b. 
As expected, it shows a wide fluctuation scattered from b = d to b = 3d. When b is 
associated with v, we find b increases upto SAd and decrease as v increases. In other 
words, we observe higher b moving away from the stagnant zone and towards the fast 
flow regions at the center. However, for the fastest-flow regions close to the orifice, 
b decreases. A reasonable implication of the increase in b is that the slightly lower 
density in the fast-flow region due to dilation makes horizontal movements easier. The 
decrease in b at higher v is perhaps related to the fact that particles undergo coUisional 
flow in the fast flowing regions near the orifice. Since the particles are less locked to 
neighbors than in the dense bulk away from the orifice, the shear in the downward 
velocity result in less horizontal movement, therefore, smaller b. 

A few further comments about figure {7\ are in order. To collect meaningful 
statistics for 6, we ignore shear-free zones (e.g. stagnant zone and plug-flow regions 
where the gradient of v is negligible), where b is likely to have large errors. We 
accomplish this by only considering cells where gradient is larger than 5% of the 
characteristic magnitude, v* /d, where v* is overall average velocity in the plug region. 
Although we only discuss b for the silo experiments in figure [T] a similar trend is also 
found for the hoppers as well. 

The correlation between v and 6 gives some clues to explain the discrepancies in 
section l^m The Kinematic Model with constant b fails to capture the development of 
a more plug-like plateau in the velocity profile even with larger values of 6. However, 
if higher b is applied to the region around the center (where v is high) , and lower b is 
applied to the region close to walls (where v is low), the model would come into closer 
agreement with experiment. 

In an effort to understand the universality of this pattern, we use the overall 
average velocity, w*, to normalize the downward velocity, v, from different flow rates (or 
orifice size). As shown in figure[7| pairs of {v/v* , b) for three different fiow rates fall into 
nearly the same pattern. This is a consistent with the trends observed in Rcf. "W that 
increasing the fiow rate merely fast-forwards the entire dynamics, without changing 
the geometrical sequence of events. 

Our way to describe our experimental results a posteriori is via a modified 
constitutive law with a variable diffusion length, b, which depends on the (scaled) 
local velocity: 



where b* is an effective diffusion length and <I> is a dimensionless scaling function. Note 
that the velocity field satisfying Eq. (|17|l is still linear with respect to rescaling the 
total magnitude of the velocity (by changing the total flow rate) since v/v* is invariant 



dx 



dv 



and 




(17) 
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when V is rescaled. However, the velocity profile in space is governed by a nonlinear 
diffusion equation, 

dv 
dx 



dv _ d 
dz dx 



/ V \ dv 
\v* ' dx 



(18) 



It is well known that spreading solutions to this equation (analogous to a 
concentration-dependent diffusivity) are flatter in the central region (compared to 
a Gaussian) when $ is an increasing function of its argument |38| . 

We should consider what might be the microscopic reason for a nonlinear diffusion 
length in the Kinematic Model. In general, it would arise from interactions between 
different spots, which are neglected as a first approximation. It makes sense that spots 
of free volume should diffuse less when they find themselves in a more slowly flowing, 
less dense, region, with fewer other nearby spots. This could explain why h appears to 
grow with velocity (or spot concentration). On the other hand, the flow in the upper 
part of the silo becomes more plug-like should exhibit less diffusion than the lower 
region of greater shear near the orifice, so it remains unclear whether the nonlinear 
model (|17|l can be given a firm microscopic justification. Further comparison with 
theory and experiment is needed to settle this question. 

The next issue to test is the boundary condition at the side walls. Specifically, it 
is important to test if the model can be simply extended from open silos to hoppers 
by using Eq. Q). It is interesting to note that the curvature of the profile aX z — zi 
around x — Q remains the same for the different funnel angles [see figure [2Ib)'-4nib)]. 
In fact, it is b that should increase from h = 1.3c? to 6 = 2.%d in order to reproduce the 
same curvature as the hopper angle is increased. For a more quantitative argument, 
we show in figure |S1 the variance of the downward velocity profile (a measure of its 
squared width). 



2, / v{x)dz 



{x'), = ■> \' (19) 
J v{x)dx 

as a function of height, z. From Eq. the slope of the linear regime near the orifice 
is equal to 26, and the value of the implied b does not significantly vary from b — 1.9c? 
for the silo, as can be seen in the inset to figure |S1 

We conclude, therefore, that extending the Kinematic Model to a hopper with 
non- vertical walls does not seem to be successful with the naive idea of Eq. Q , which 
assumes the same bulk constitutive law holds all the way to the boundary. 

It may be that a nonlinear constitutive law as in Eq. H17|l can improve the situation 
because particles slip more on a funnel wall and b thus tends to be higher than in the 
silo. However, there may still be problems higher in the tank where the fiowing region 
meets the vertical side walls. The boundary condition requires that the strain 
rate (horizontal gradient of vertical velocity) vanishes at a vertical wall, and yet small 
velocity gradients are observed near the walls in the upper region in Figs. 2-5. We 
plan to compare the nonlinear Kinematic Model, as well as other continuum models 
from critical-state mechanics and hourglass theory, more closely with the experimental 
flow profiles in future work. 



6. Summary 

In summary, we have used high speed imaging techniques to track the positions of 
granular materials draining inside silos and hoppers. We compared our data with the 
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Figure 8. (color online). The variance (squared effective width) of the downward 
velocity profile, (x'^)v, as a function of the vertical coordinate, z for different funnel 
angles. 



continuum Kinetic Model and two possible microscopic theories which predict similar 
mean flow, the Void Model and Spot Model. These models are appealing due to 
their mathematical simplicity and completeness, which allows direct application to 
various geometries. The models also predict smooth velocity profiles, free of shock- 
like discontinuities, quite consistent with the experiments. Systematic deviations are 
observed, implying various assumptions, such as a constant diffusion length, are too 
simple to capture all aspects of the flow profile, but it may be that modifications can be 
made to improve the agreement. For example, we infer that the kinematic parameter, 
6, increases with the local velocity, which would imply that the spot diffusion length 
increases in the presence of other spots. Still, it is clear that further work is also 
needed to develop boundary conditions for both discrete and continuous models of 
slow, dense granular flows. 
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